function [stable_teach_propose,stable_school_propose,firstbest_school,teach_school_same,jass_iprop_char, jass_jprop_char] = func_allocations( IJipref,IJipref_choice,IJoutput, IJVAma,IJclassSize,IJVAm1, IJVAm2,IJVAmean, IJfracDisadv, IJitidx,IJjtidx,IJblack,IJexp7p,IJappyear,yearselect,makebalanced,errors_on,firstbest)

% This takes in information about teacher and school preferences.
% Solves for (1) teacher propose stable matches; (2) school propose stable
% matching (3) allocations that maximize school objectives

% In the two school objectives, computes characteristics of the teachers to
% which each school is assigned

% equilibrium allocation output:
%stable_teach_propose=[stable_teach_propose_util stable_teach_propose_achiev stable_teach_propose_type1 stable_teach_propose_type2];
%stable_school_propose=[stable_school_propose_util stable_school_propose_achiev stable_school_propose_type1 stable_school_propose_type2];
%firstbest_school=[firstbest_school_util firstbest_school_achiev firstbest_school_type1 firstbest_school_type2];

% characteristics of teacher to which school is assigned (i is teacher, and j i school)
%jass_iprop_char=[IJfracDisadv_ij(:,1) IJVAmean_sq(jass_iprop,1) IJblack_sq(jass_iprop,1)  IJexp7p_sq(jass_iprop,1)  IJVAca_sq(jass_iprop,1) ]; 
%jass_jprop_char=[IJfracDisadv_ij(:,1) IJVAmean_sq(jass_jprop,1) IJblack_sq(jass_jprop,1)  IJexp7p_sq(jass_jprop,1)  IJVAca_sq(jass_jprop,1) ]; 
  


% User supplies the year of the data to use
% User specifies whether make the schools and teachers balanced
    %% go down to IJappyear==yearselect

    
    IJipref=IJipref(IJappyear==yearselect);
    IJipref_choice=IJipref_choice(IJappyear==yearselect);
    IJitidx2015=IJitidx(IJappyear==yearselect);
    IJjtidx2015=IJjtidx(IJappyear==yearselect);
    IJVAma=IJVAma(IJappyear==yearselect);
    IJclassSize=IJclassSize(IJappyear==yearselect);
    IJoutput=IJoutput(IJappyear==yearselect);
    
    
    IJVAm1=IJVAm1(IJappyear==yearselect);
    IJVAm2=IJVAm2(IJappyear==yearselect);
    IJfracDisadv=IJfracDisadv(IJappyear==yearselect);
    IJVAmean=IJVAmean(IJappyear==yearselect);
    
    IJblack=IJblack(IJappyear==yearselect);
    IJexp7p=IJexp7p(IJappyear==yearselect);
    
    % redo indexes so that they start at 1 
    
    IJitidx=IJitidx2015-min(IJitidx2015)+1;
    IJjtidx=IJjtidx2015-min(IJjtidx2015)+1;
    
    if makebalanced==1
        imax=max(IJitidx);       
        jlist=unique(IJjtidx);
        jkeep=randperm(max(jlist),imax)';  % select a number of classrooms that equals the number of teachers
        jselect=ismember(IJjtidx,jkeep);  % generate a logical vector telling us which to keep
        
        % now apply the selector
        
        IJipref=IJipref( jselect);
        IJipref_choice=IJipref_choice( jselect);       
        IJVAma=IJVAma( jselect);
        IJclassSize=IJclassSize( jselect);
        IJoutput=IJoutput( jselect);
        IJVAmean=IJVAmean(jselect);
        IJVAm1=IJVAm1( jselect);
        IJVAm2=IJVAm2( jselect);
        IJfracDisadv= IJfracDisadv( jselect);
        IJitidx=IJitidx( jselect);
        IJjtidx= IJjtidx( jselect);
        IJblack=IJblack(jselect);
        IJexp7p=IJexp7p(jselect);       
        
        % sort j in sequential order (to renumber), then reorder everything
        % else
        [b1,Isorter]=sort(IJjtidx); % sort by j
        inbetween=repmat((1:imax),imax,1);
        IJjtidx=inbetween(:);
        
        IJipref=IJipref(Isorter);
        IJipref_choice=IJipref_choice(Isorter);        
        IJVAma=IJVAma( Isorter);
        IJclassSize=IJclassSize( Isorter);
        IJoutput=IJoutput( Isorter);
        IJVAmean=IJVAmean(Isorter);
        IJVAm1=IJVAm1( Isorter);
        IJVAm2=IJVAm2(Isorter);
        IJfracDisadv= IJfracDisadv(Isorter);
        IJitidx=IJitidx(Isorter);
        IJblack=IJblack(Isorter);
        IJexp7p=IJexp7p(Isorter);      

    end
%% make the teacher preferences rectangular, and shocks, and then rank
    
    error_draws=evrnd(0,1,length(IJitidx),1);
    
    if errors_on==1
            for i=1:length(IJitidx)
                 IJipref_sq(IJitidx(i,1),IJjtidx(i,1))=IJipref(i,1)+error_draws(i,1); % add the type I errors
                 IJipref_choice_sq(IJitidx(i,1),IJjtidx(i,1))=IJipref_choice(i,1)+error_draws(i,1); % add the type I errors
            end
    end
    
    if errors_on==0
         for i=1:length(IJitidx)
            IJipref_sq(IJitidx(i,1),IJjtidx(i,1))=IJipref(i,1)+error_draws(i,1); % even if errors aren't on, need to make the experienced utility have the errors
            IJipref_choice_sq(IJitidx(i,1),IJjtidx(i,1))=IJipref_choice(i,1);   % even if errors aren't on, need to make this square          
         end
    end

    [B1, IJirank]=sort(IJipref_choice_sq,2,'descend');  % this is sorted in descending order; so the first column in I is the most preferred of the teacher in the relevant row 
    [B1, IJirank_exp]=sort(IJipref_sq,2,'descend');  % this is sorted in descending order; so the first column in I is the most preferred of the teacher in the relevant row 

    % make the school output vector rectangular, and then rank
    % note: this breaks ties in an arbitrary way 
   
   for i=1:length(IJitidx)
         IJjpref_sq(IJjtidx(i,1),IJitidx(i,1))=IJoutput(i,1);    % these are J X I (useful for DA)
         IJVAma_sq(IJjtidx(i,1),IJitidx(i,1))=IJVAma(i,1); % these are J X I (useful for DA)
         IJfracDisadv_ij(IJjtidx(i,1),IJitidx(i,1))=IJfracDisadv(i,1);  % these are J x 1 (useful for first best)   
         IJclassSize_ij(IJjtidx(i,1),IJitidx(i,1))=IJclassSize(i,1);  % these are I x J (useful for first best)
         IJijpref_sq(IJitidx(i,1),IJjtidx(i,1))=IJoutput(i,1);    % these are I X J (useful for first best)
         IJijVAma_sq(IJitidx(i,1),IJjtidx(i,1))=IJVAma(i,1);  % these are I X J (useful for first best)  
         IJclassSize_sq(IJitidx(i,1),IJjtidx(i,1))=IJclassSize(i,1);  % these are I x J (useful for first best)
         IJVAm1_sq(IJitidx(i,1),IJjtidx(i,1))=IJVAm1(i,1);  % these are I X J (useful for first best)
         IJVAm2_sq(IJitidx(i,1),IJjtidx(i,1))=IJVAm2(i,1);  % these are I X J (useful for first best)
         IJVAmean_sq(IJitidx(i,1),IJjtidx(i,1))=IJVAmean(i,1); % I x J 
         IJfracDisadv_sq(IJitidx(i,1),IJjtidx(i,1))=IJfracDisadv(i,1);  % these are I X J (useful for first best)    
         IJblack_sq(IJitidx(i,1),IJjtidx(i,1))=IJblack(i,1);
         IJexp7p_sq(IJitidx(i,1),IJjtidx(i,1))=IJexp7p(i,1);      
   end
   
   IJtype1num=(1-IJfracDisadv_sq).*IJclassSize_sq;
   IJtype2num=IJfracDisadv_sq.*IJclassSize_sq;
       
   [B1, IJjrank]=sort(IJjpref_sq,2,'descend');  % this is sorted in descending order; so the first column in I is the most preferred of the school in the relevant row 
    
   % run DA
    
   [iass_jprop]= galeshapley_gender_unbalanced(IJjrank,IJirank);  % school propose
   
   %  figure out jass_jprop (including the 0s)
   
   jstuff= [ (1:length(iass_jprop))' iass_jprop];
   jstuff2=sortrows(jstuff,2);
   
   jmax=max(IJjtidx);
   imax=max(IJitidx);
   
   jass_jprop=[ zeros(jmax,1)];
   for jj=1:jmax
       for j=1:imax
         if (jstuff2(j,2)==jj)
            jass_jprop(jj,1)=jstuff2(j,1);
         end
       end
   end
          
   [jass_iprop]= galeshapley_gender_unbalanced(IJirank,IJjrank);  % teach propose
    
   % figure out iass_iprop
   istuff= [ (1:length(jass_iprop))' jass_iprop];
   istuff2=sortrows(istuff,2);
      
   iass_iprop=[ zeros(imax,1)];
   for ii=1:imax
       for i=1:jmax
         if (istuff2(i,2)==ii)
            iass_iprop(ii,1)=istuff2(i,1);
         end
       end
   end     
   
   
  % record characteristics of the teachers to which the school is assigned,
  % and one school char
  IJVAca_sq=IJVAm2_sq-IJVAm1_sq;
  jass_iprop_char=[IJfracDisadv_ij(:,1) IJVAmean_sq(jass_iprop,1) IJblack_sq(jass_iprop,1) IJexp7p_sq(jass_iprop,1)  IJVAca_sq(jass_iprop,1) ]; 
  jass_jprop_char=[IJfracDisadv_ij(:,1) IJVAmean_sq(jass_jprop,1) IJblack_sq(jass_jprop,1) IJexp7p_sq(jass_jprop,1)  IJVAca_sq(jass_jprop,1) ]; 
  
  % test if the school propose and teacher propose are the same
  
  teach_school_same=min((iass_iprop==iass_jprop));

  % figure out if a teacher is unassigned, if so, remove the row from the matrices
  % expand this to do it separately for teach and school propose 
   keep_teachj=logical(iass_jprop>0);
   IJipref_sqj=IJipref_sq(keep_teachj,:);
   IJVAm1_sqj=IJVAm1_sq(keep_teachj,:);
   IJtype1numj=IJtype1num(keep_teachj,:);
   IJVAm2_sqj= IJVAm2_sq(keep_teachj,:);
   IJtype2numj=IJtype2num(keep_teachj,:);
   IJclassSize_sqj=IJclassSize_sq(keep_teachj,:);
   IJVAmean_sqj=IJVAmean_sq(keep_teachj,:);
   iass_jprop=iass_jprop(keep_teachj,:);
   IJijVAma_sqj=IJijVAma_sq(keep_teachj,:);
   IJirankj=IJirank(keep_teachj,:);
   IJirank_expj=IJirank_exp(keep_teachj,:);

   IJfracDisadv_sqj=IJfracDisadv_sq(keep_teachj,:);
   
   keep_teachi=logical(iass_iprop>0);
   IJipref_sqi=IJipref_sq(keep_teachi,:);
   IJVAm1_sqi=IJVAm1_sq(keep_teachi,:);
   IJtype1numi=IJtype1num(keep_teachi,:);
   IJVAm2_sqi= IJVAm2_sq(keep_teachi,:);
   IJtype2numi=IJtype2num(keep_teachi,:);
   IJclassSize_sqi=IJclassSize_sq(keep_teachi,:);
   IJVAmean_sqi=IJVAmean_sq(keep_teachi,:);
   iass_iprop=iass_iprop(keep_teachi,:);
   IJijVAma_sqi=IJijVAma_sq(keep_teachi,:);
   IJiranki=IJirank(keep_teachi,:);
   IJirank_expi=IJirank_exp(keep_teachi,:);
   IJfracDisadv_sqi=IJfracDisadv_sq(keep_teachi,:);
   
   % compute per teacher utility in teacher propose and school propose

   for i=1:min(max(IJitidx),max(IJjtidx))
     i_util_iprop(i,1)=IJipref_sqi(i,iass_iprop(i,1));
   end
   
   for i=1:min(max(IJitidx),max(IJjtidx))
     i_util_jprop(i,1)=IJipref_sqj(i,iass_jprop(i,1));
   end
    
   % compute type I average output, and type 2 average output in teach
   % propose
   for i=1:min(max(IJitidx),max(IJjtidx))
          i_type1_tot(i,1)=IJVAm1_sqi(i,iass_iprop(i,1)).* IJtype1numi(i,iass_iprop(i,1));
          i_type1_stu(i,1)= IJtype1numi(i,iass_iprop(i,1));
          
          i_type2_tot(i,1)=IJVAm2_sqi(i,iass_iprop(i,1)).* IJtype2numi(i,iass_iprop(i,1));
          i_type2_stu(i,1)= IJtype2numi(i,iass_iprop(i,1));
   end
        
   stable_teach_propose_type1=sum(i_type1_tot)/sum(i_type1_stu);
   stable_teach_propose_type2=sum(i_type2_tot)/sum(i_type2_stu);
       
   % compute type I average output, and type 2 average output in
   % school propose
   for i=1:min(max(IJitidx),max(IJjtidx))
          i_type1_tot(i,1)=IJVAm1_sqj(i,iass_jprop(i,1)).* IJtype1numj(i,iass_jprop(i,1));
          i_type1_stu(i,1)= IJtype1numj(i,iass_jprop(i,1));
          
          i_type2_tot(i,1)=IJVAm2_sqj(i,iass_jprop(i,1)).* IJtype2numj(i,iass_jprop(i,1));
          i_type2_stu(i,1)= IJtype2numj(i,iass_jprop(i,1));
   end
        
   stable_school_propose_type1=sum(i_type1_tot)/sum(i_type1_stu);
   stable_school_propose_type2=sum(i_type2_tot)/sum(i_type2_stu);  
   type2_stu_iprop=i_type2_stu;
  
   
   % compute per student output in each assignment 
   % loop over teachers: all the teachers are assigned whereas the classes aren't
    
  for i=1:min(max(IJitidx),max(IJjtidx))
          i_output_iprop(i,1)=IJijVAma_sqi(i,iass_iprop(i,1)).*IJclassSize_sqi(i,iass_iprop(i,1));
          i_num_students(i,1)=IJclassSize_sqi(i,iass_iprop(i,1));

  end
  
  stable_teach_propose_achiev=sum(i_output_iprop)/sum(i_num_students);
   
  for i=1:min(max(IJitidx),max(IJjtidx))
         i_output_jprop(i,1)=IJijVAma_sqj(i,iass_jprop(i,1)).*IJclassSize_sqj(i,iass_jprop(i,1));
         i_num_students(i,1)=IJclassSize_sqj(i,iass_jprop(i,1));
  end
    
  stable_school_propose_achiev=sum(i_output_jprop)/sum(i_num_students);

  % compute first best 
  % restrict to the classrooms that are filled (select columns) 
  % use the i-j organized things for output 
  
  if firstbest==1 
  % student achievement maximizing (or: school objective maximizing)
  out_sq= -IJijpref_sq(keep_teachi,iass_iprop);
  va_output_sq=IJijpref_sq(keep_teachi,iass_iprop);
  va_sq = IJijVAma_sq(keep_teachi,iass_iprop);
  teachpref_sq=IJipref_choice_sq(keep_teachi,iass_iprop);
  teachpref_true_sq=IJipref_sq(keep_teachi,iass_iprop);
  IJijVAma_sq=IJijVAma_sq(keep_teachi,iass_iprop);
  IJVAm1_sq=IJVAm1_sq(keep_teachi,iass_iprop);
  IJVAm2_sq=IJVAm2_sq(keep_teachi,iass_iprop);
  IJtype1num=IJtype1num(keep_teachi,iass_iprop); 
  IJtype2num=IJtype2num(keep_teachi,iass_iprop); 
  IJclassSize_sq=IJclassSize_sq(keep_teachi,iass_iprop);  
  IJfracDisadv_sq=IJfracDisadv_sq(keep_teachi,iass_iprop);
  IJVAmean_sq=IJVAmean_sq(keep_teachi,iass_iprop);


  
  IJjrank=IJjrank(iass_iprop,:);
    
  [iind_best,jind_best] = linear_sum_assignment(out_sq);

  student_max_achievement=zeros(length(iind_best),1);
  student_max_pref=zeros(length(iind_best),1);
for i=1:length(iind_best);
    student_max_pref(i,1)= teachpref_sq(iind_best(i,1),jind_best(i,1));
end

   % compute type I average output, and type 2 average output in max school
        for i=1:min(max(IJitidx),max(IJjtidx))      
          i_type1_tot(i,1)=IJVAm1_sq(iind_best(i,1),jind_best(i,1)).* IJtype1num(iind_best(i,1),jind_best(i,1));
          i_type1_stu(i,1)= IJtype1num(iind_best(i,1),jind_best(i,1));
          
          i_type2_tot(i,1)=IJVAm2_sq(iind_best(i,1),jind_best(i,1)).* IJtype2num(iind_best(i,1),jind_best(i,1));
          i_type2_stu(i,1)= IJtype2num(iind_best(i,1),jind_best(i,1));
          
          i_over_tot(i,1)=  IJijVAma_sq(iind_best(i,1),jind_best(i,1)).*IJclassSize_sq(iind_best(i,1),jind_best(i,1));
          i_over_stu(i,1)=IJclassSize_sq(iind_best(i,1),jind_best(i,1));     
        end
        
        firstbest_school_type1=sum(i_type1_tot)/sum(i_type1_stu);
        firstbest_school_type2=sum(i_type2_tot)/sum(i_type2_stu);
        firstbest_school_achiev=sum(i_over_tot)/sum(i_over_stu);

 firstbest_school_util=mean(student_max_pref);
firstbest_school=[firstbest_school_util firstbest_school_achiev firstbest_school_type1 firstbest_school_type2];       
  end 

if firstbest==0
    firstbest_school=zeros(1,4);
end
  
 
% collect output

stable_teach_propose_util= mean(i_util_iprop);
stable_school_propose_util= mean(i_util_jprop);


stable_teach_propose=[stable_teach_propose_util stable_teach_propose_achiev stable_teach_propose_type1 stable_teach_propose_type2];
stable_school_propose=[stable_school_propose_util stable_school_propose_achiev stable_school_propose_type1 stable_school_propose_type2];


end

